COMMENT ⊗   VALID 00006 PAGES
C REC  PAGE   DESCRIPTION
C00001 00001
C00002 00002	\input acphdr % Section 4.7
C00003 00003	%folio 635 galley 5b (C) Addison-Wesley 1978	*
C00012 00004	%folio 638 galley 6 (C) Addison-Wesley 1978	*
C00021 00005	%folio 643 galley 7 (C) Addison-Wesley 1978	*
C00044 00006	\eject
C00045 ENDMK
C⊗;
\input acphdr % Section 4.7
\open0=v247.inx
\runninglefthead{ARITHMETIC} % chapter title
\titlepage\setcount00
\null
\vfill
\tenpoint
\ctrline{SECTION 4.7 of THE ART OF COMPUTER PROGRAMMING}
\ctrline{$\copyright$ 1980 Addison--Wesley Publishing Company, Inc.}
\vfill
\runningrighthead{MANIPULATION OF POWER SERIES}
\section{4.7}
\eject
\setcount0 486
%folio 635 galley 5b (C) Addison-Wesley 1978	*
\sectionbegin{\star4.7. MANIPULATION OF POWER SERIES}
I{\:mF WE ARE GIVEN} two power series
$$U(z) = U↓0 + U↓1z + U↓2z↑2 +\cdotss,\qquad V(z) = V↓0 + V↓1z
+ V↓2z↑2 +\cdotss,\eqno (1)$$
whose coefficients belong to a field, we can form
their sum, their product, their quotient, etc., to obtain new
power series. A polynomial is obviously a special case of a
power series, in which there are only finitely many terms.

Of course, only a finite number of terms can be
represented and stored within a computer, so it makes sense
to ask whether power series arithmetic is even possible on computers;
and if it is possible, what makes it different from polynomial
arithmetic? The answer is that we work with only the first $N$
coefficients of the power series, where $N$ is a parameter that
may in principle be arbitrarily large; instead of ordinary polynomial
arithmetic, we are essentially doing polynomial arithmetic modulo
$z↑N$, and this often leads to a seomewhat different point of
view. Furthermore, special operations like ``reversion'' can
be performed on power series but not on polynomials, since polynomials
are not closed under these operations.

Manipulation of power series has several applications
to numerical analysis, but perhaps its greatest use is the determination
of asymptotic expansions (as we have seen in Section 1.2.11.3),
or the calculation of quantities defined by certain generating
functions. The latter applications make it desirable to calculate
the coefficients exactly, instead of with floating-point arithmetic.
All of the algorithms in this section, with obvious exceptions,
can be done using rational operations only, so the techniques
of Section 4.5.1 can be used to obtain exact results when desired.

The calculation of $W(z) = U(z) \pm V(z)$ is, of
course, trivial, since we have $W↓n = U↓n \pm V↓n$ for $n =
0$, 1, 2, $\ldotss\,$. It is also easy to calculate $W(z) = U(z)V(z)$,
using the familiar ``Cauchy product rule'':
$$\chop to 9pt{W↓n = \sum ↓{0≤k≤n} U↓kV↓{n-k} = U↓0V↓n + U↓1V↓{n-1} +\cdots
+ U↓nV↓0.}\eqno (2)$$

The quotient $W(z) = U(z)/V(z)$, when $V↓0
≠ 0$, can be obtained by interchanging $U$ and $W$ in (2); we
obtain the rule
$$\eqalignno{W↓n ⊗= \bigglp U↓n - \sum ↓{0≤k<n} W↓kV↓{n-k}\biggrp\vcenter{
\hbox{\:@\char'36}}V↓0\cr
⊗= (U↓n - W↓0V↓n - W↓1V↓{n-1} -\cdots - W↓{n-1}V↓1)/V↓0.⊗(3)\cr}$$
This recurrence relation for the $W$'s makes it easy
to determine $W↓0$, $W↓1$, $W↓2$, $\ldots$ successively, without inputting
$U↓n$ and $V↓n$ until after $W↓{n-1}$ has been computed. Let
us say that a power-series manipulation algorithm with the latter
property is ``on-line''; an on-line algorithm can be used to
determine $N$ coefficients $W↓0$, $W↓1$, $\ldotss$, $W↓{N-1}$ of the
result without knowing $N$ in advance, so it is possible in
theory to run the algorithm indefinitely and compute the entire
power series; or to run it until a certain condition is met.\xskip
(The opposite of ``on-line'' is ``off-line.'')

If the coefficients $U↓k$ and $V↓k$ are integers
but the $W↓k$ are not, recurrence (3) involves computation with
fractions. This can be avoided by the all-integer approach described
in exercise 2.

Let us now consider the operation of computing
$W(z) = V(z)↑α$, where $α$ is an ``arbitrary'' power. For
example, we could calculate the square root of $V(z)$ by taking
$α = {1\over 2}$, or we could find $V(z)↑{-10}$ or even $V(z)↑π$.
If $V↓m$ is the first nonzero coefficient of $V(z)$, we have
$$\baselineskip15pt\eqalign{V(z)⊗= V↓{m\,}z↑m\biglp 1 + (V↓{m+1}/V↓m)z
+ (V↓{m+2}/V↓m)z↑2 +\cdotss\bigrp ,\cr
V(z)↑α ⊗= V↑{α}↓{\!m\,}z↑{αm}\biglp 1 + (V↓{m+1}/V↓m)z + (V↓{m+2}/V↓m)z↑2
+\cdotss\bigrp ↑α.\cr}\eqno(4)$$
This will be a power series if and only if $αm$ is a
nonnegative integer. From (4) we can see that the problem of
computing general powers can be reduced to the case that $V↓0
= 1$; then the problem is to find coefficients of
$$W(z) = (1 + V↓1z + V↓2z↑2 + V↓3z↑3 +\cdotss)↑α.\eqno(5)$$
Clearly $W↓0 = 1↑α = 1$.

The obvious way to find the coefficients of (5)
is to use the binomial theorem (Eq.\ 1.2.9--19), or (if $α$ is
a positive integer) to try repeated squaring as in Section 4.6.3;
but a much simpler and more efficient device for calculating
powers has been suggested by J. C. P. Miller.\xskip [See P. Henrici,
{\sl JACM \bf 3} (1956), 10--15.]\xskip If $W(z) = V(z)↑α$, we have
by differentiation
$$W↓1 + 2W↓2z + 3W↓3z↑2 +\cdots = W↑\prime (z) =
αV(z)↑{α-1}V↑\prime(z);\eqno (6)$$
therefore
$$W↑\prime(z)V(z) = αW(z)V↑\prime(z).\eqno (7)$$
If we now equate the coefficients of $z↑{n-1}$
in (7), we find that
$$\sum ↓{0≤k≤n}kW↓kV↓{n-k} = α\sum ↓{0≤k≤n} (n - k)W↓kV↓{n-k},\eqno
(8)$$
and this gives us a useful computational rule valid for all $n≥1$:
$$\eqalignno{W↓n ⊗= \sum ↓{1≤k≤n} \left(\left(α + 1\over n\right)\,k
- 1\right)\,V↓kW↓{n-k}\cr
⊗= \biglp (α{+}1{-}n)V↓1W↓{n-1} + (2α{+}2{-}n)V↓2W↓{n-2}
+\cdots + nαV↓n\bigrp /n.⊗(9)\cr}$$
This equation leads to a simple on-line algorithm by which we can successively
determine $W↓1$, $W↓2$, $\ldotss$, using approximately $2n$ multiplications to
compute the $n$th coefficient. Note the special case $α=-1$, in which (9)
becomes the special case $U(z)=V↓0=1$ of (3).

A similar technique can be used to form $f\biglp V(z)\bigrp$ when $f$ is any
function that satisfies a simple differential equation.\xskip(For example, see
exercise 4.)\xskip A comparatively straightforward ``power-series method'' is
often used to obtain the solution of differential equations; this technique is
explained in nearly all textbooks about differential equations.
%folio 638 galley 6 (C) Addison-Wesley 1978	*
\subsectionbegin{Reversion of series}
The transformation of power series that
is perhaps of greatest interest is called ``reversion of series.''
This problem is to solve the equation
$$z = t + V↓2t↑2 + V↓3t↑3 + V↓4t↑4 +\cdots \eqno(10)$$
for $t$, obtaining the coefficients of the power series
$$t = z + W↓2z↑2 + W↓3z↑3 + W↓4z↑4 +\cdotss.\eqno (11)$$

Several interesting ways to achieve such a reversion
are known. We might say that the ``classical'' method is one
based on Lagrange's remarkable inversion formula [{\sl M\'emoires
Acad.\ Royale des Sciences et Belles-Lettres de Berlin \bf 24}
(1768), 251--326], which states that
$$W↓n = U↓{n-1}/n,$$
if
$$U↓0 + U↓1t + U↓2t↑2 +\cdots = (1 + V↓2t +
V↓3t↑2 +\cdotss)↑{-n}.\eqno (12)$$
For example, we have $(1 - t)↑{-5} = {4\choose 4}
+ {5\choose 4}t + {6\choose 4}t↑2 +\cdotss$; hence $W↓5$ in
the reversion of $z = t - t↑2$ is equal to ${8\choose 4}/5
= 14$. This checks with the formulas for enumerating binary trees
in Section 2.3.4.4.

Relation (12) shows that we can revert the series
(10) if we compute the negative powers $(1 + V↓2t + V↓3t↑2 +\cdotss)↑{-n}$
for $n = 1$, 2, 3, $\ldotss\,$. A straightforward application of
this idea would lead to an on-line reversion algorithm that
uses approximately $N↑3/2$ multiplications to find $N$ coefficients,
but Eq.\ (9) makes it possible to work with only the first $n$
coefficients of $(1 + V↓2t + V↓3t↑2 +\cdotss)↑{-n}$,
obtaining an on-line algorithm that requires only about $N↑3/6$
multiplications.

\algbegin Algorithm L (Lagrangian power-series
reversion). This on-line algorithm inputs the value of $V↓n$
in (10) and outputs the value of $W↓n$ in (11), for $n = 2$,
3, 4, $\ldotss$,\penalty999\ $N$.\xskip (The number $N$ need not be specified in
advance; some other termination condition may be substituted.)

\algstep L1. [Initialize.] Set $n ← 1$, $U↓0
← 1$.\xskip (The relation
$$(1 + V↓2t + V↓3t↑2 +\cdotss)↑{-n}
= U↓0 + U↓1t +\cdots + U↓{n-1}t↑{n-1} + O(t↑n)\eqno(13)$$
will be maintained throughout this algorithm.)

\topinsert{\vskip 34mm
\ctrline{\caption Fig.\ 16. Power-series reversion by Algorithm L.}}

\algstep L2. [Input $V↓n$.] Increase $n$ by 1. If
$n > N$, the algorithm terminates; otherwise input the next
coefficient, $V↓n$.

\algstep L3. [Divide.] Set $U↓k ← U↓k - U↓{k-1}V↓2 -\cdots
- U↓1V↓k - V↓{k+1}$, for $k = 1$, 2, $\ldotss$, $n - 2$ (in
this order); then set
$$U↓{n-1} ← -2U↓{n-2}V↓2 - 3U↓{n-3}V↓3 -\cdots
- (n - 1)U↓1V↓{n-1} - nV↓{n-1}.$$
$\biglp$We have thereby divided $U(z)$ by
$V(z)/z$; cf.\ (3) and (9).$\bigrp$

\algstep L4. [Output $W↓n$.] Output $U↓{n-1}/n$
(which is $W↓n$) and return to L2.\quad\blackslug

\yyskip When applied to the example $z = t - t↑2$,
the computation in Algorithm L is
$$\vbox{\halign{$\ctr{#}$\qquad⊗$\hfill#$\qquad⊗$\hfill#$\qquad⊗$\hfill#$\qquad
⊗$\hfill#$\qquad⊗$\hfill#$\qquad⊗$\hfill#$\qquad⊗$\hfill#$\cr
n⊗V↓n\hfill⊗U↓0\hfill⊗U↓1\hfill⊗U↓2\hfill⊗U↓3\hfill⊗U↓4\hfill⊗W↓n\cr
\noalign{\vskip3pt}
1⊗1⊗1⊗⊗⊗⊗⊗1\cr
2⊗-1⊗1⊗2⊗⊗⊗⊗1\cr
3⊗0⊗1⊗3⊗6⊗⊗⊗2\cr
4⊗0⊗1⊗4⊗10⊗20⊗⊗5\cr
5⊗0⊗1⊗5⊗15⊗35⊗70⊗14\cr}}$$
Exercise 8 shows that a slight modification of Algorithm
L will solve a considerably more general problem with only a
little more effort.

\yskip Let us now consider solving the equation
$$U↓1z + U↓2z↑2 + U↓3z↑3 +\cdots = t + V↓2t↑2 +
V↓3t↑3 +\cdots \eqno (14)$$
for $t$, obtaining the coefficients of the power
series
$$t = W↓1z + W↓2z↑2 + W↓3z↑3 + W↓4z↑4 +\cdotss.\eqno (15)$$
Eq.\ (10) is the special case $U↓1 = 1$, $U↓2 = U↓3
=\cdots = 0$. If $U↓1 ≠ 0$, we may assume that $U↓1
= 1$, if we replace $z$ by $(U↓1z)$; but we shall consider the
general equation (14), since $U↓1$ might equal zero.

\algbegin Algorithm T (General power-series reversion).
This on-line algorithm inputs the values of $U↓n$ and $V↓n$
in (14) and outputs the value of $W↓n$ in (15), for $n = 1$,
2, 3, $\ldotss$, $N$. An auxiliary matrix $T↓{mn}$, $1 ≤ m ≤ n ≤
N$, is used in the calculations.

\algstep T1. [Initialize.] Set $n ← 1$. Let
the first two inputs (namely, $U↓1$ and $V↓1$) be stored in
$T↓{11}$ and $V↓1$, respectively.\xskip (We must have $V↓1 = 1.$)

\algstep T2. [Output $W↓n$.] Output the value of $T↓{1n}$ (which
is $W↓n$).

\algstep T3. [Input $U↓n$, $V↓n$.] Increase $n$ by 1. If $n
> N$, the algorithm terminates; otherwise store the next two
inputs (namely $U↓n$ and $V↓n$) in $T↓{1n}$ and $V↓n$, respectively.

\algstep T4. [Multiply.] Set
$$T↓{mn} ← T↓{11}T↓{m-1,n-1} + T↓{12}T↓{m-1,n-2} +\cdots
+ T↓{1,n-m+1}T↓{m-1,m-1}$$
and $T↓{1n} ← T↓{1n} - V↓mT↓{mn}$,
for $2 ≤ m ≤ n$.\xskip $\biglp$After this step we have the conditions
$$t↑m = T↓{mm}z↑m + T↓{m,m+1}z↑{m+1} +\cdots + T↓{mn}z↑n
+ O(z↑{n+1}),\qquad 1 ≤ m ≤ n.\eqno(16)$$
It is easy to verify (16) by induction
for $m ≥ 2$, and when $m = 1$, we have $U↓n = T↓{1n} + V↓2T↓{2n}
+\cdots + V↓nT↓{nn}$ by (14) and (19).$\bigrp$\xskip Return to
step T2.\quad\blackslug

\yyskip Equation (16) explains the
mechanism of this algorithm, which is due to Henry C. Thacher,
Jr.\ [{\sl CACM \bf 9} (1966), 10--11]. The running time is essentially
the same as Algorithm L\null, but considerably more storage space
is required. An example of this algorithm is worked in exercise 9.
%folio 643 galley 7 (C) Addison-Wesley 1978	*
\yskip Still another approach to power series
reversion has been proposed by R.\penalty999\ P. Brent and H. T. Kung [{\sl
JACM \bf25} (1978), 581--595], who observed that standard
iterative procedures used to find roots of equations over the
real numbers can also be applied to equations over
power series. In particular, we can consider Newton's method for computing
approximations to a real number $t$ such that $f(t) = 0$, given
a function $f$ that is well-behaved near $t$: If $x$ is a good
approximation to $T$, then $\phi (x) = x - f(x)/f↑\prime (x)$
will be even better, for if we write $x = t + ε$ we have
$f(x) = f(t) + ε f↑\prime (t) + O(ε ↑2)$, $f↑\prime (x)
= f↑\prime (t) + O(ε )$, and $\varphi (x) = t + ε -
\biglp 0 + ε f↑\prime (t) + O(ε ↑2)\bigrp /\biglp f↑\prime
(t) + O(ε )\bigrp = t + O(ε ↑2)$. Applying this idea
to power series, let $f(x) = V(x) - U(z)$, where $U$ and $V$
are the power series in Eq.\ (14). We wish to find the power series
$t$ in $z$ such that $f(t) = 0$. Let $x = W↓1z +\cdots + W↓{n-1}z↑{n-1}
= t + O(z↑n)$ be an ``approximation'' to $t$ of order $n$; then
$\varphi (x) = x - f(x)/f↑\prime (x)$ will be an approximation
of order $2n$, since the assumptions of Newton's method hold
for this $f$ and $t$.

In other words, we can use the following procedure:

\algbegin Algorithm N (General power series reversion by Newton's
method). This ``semi-on-line'' algorithm inputs the
values of $U↓n$ and $V↓n$ in (14) for $2↑k ≤ n < 2↑{k+1}$ and
then outputs the values of $W↓n$ in (15) for $2↑k ≤ n < 2↑{k+1}$,
thereby producing its answers in batches of $2↑k$ at a time,
for $k = 0$, 1, 2, $\ldotss$, $K$.

\algstep N1. [Initialize.] Set $N ← 1$.\xskip (We
will have $N = 2↑k$.)\xskip Input the first coefficients $U↓1$ and
$V↓1$ (where $V↓1 = 1$), and set $W↓1 ← U↓1$.

\algstep N2. [Output.] Output $W↓n$ for $N ≤ n < 2N$.

\algstep N3. [Input.] Set $N ← 2N$. If $N > 2↑K$, the algorithm
terminates; otherwise input the values $U↓n$ and $V↓n$ for $N
≤ n < 2N$.

\algstep N4. [Newtonian step.] Use an algorithm for power series
composition (see exercise 11) to evaluate the coefficients $Q↓j$ and
$R↓j$ ($0 ≤ j < N$) in the power series
$$\cpile{\hbox to 329pt{$\dispstyle U↓1z +\cdots + U↓{2N-1}z↑{2N-1} - V(W↓1z
+\cdots + W↓{N-1}z↑{N-1})$\hfill}\cr
\noalign{\vskip2pt}
\hbox to 329pt{\hfill$\dispstyle= R↓0z↑N 
+ R↓1z↑{N+1} +\cdots + R↓{N-1}z↑{2N-1}+ O(z↑{2N}),$}\cr
\noalign{\vskip4pt}
V↑\prime (W↓1z +\cdots + W↓{N-1}z↑{N-1})=
Q↓0 + Q↓1z +\cdots + Q↓{N-1}z↑{N-1} + O(z↑N),\cr}$$
where $V(x) = x + V↓2x↑2 +\cdots$
and $V↑\prime (x) = 1 + 2V↓2x +\cdotss$. Then set $W↓N$, $\ldotss$,
$W↓{2N-1}$ to the coefficients in the power series
$$ {R↓0{+}R↓1z{+}\cdots{+}R↓{N-1}z↑{N-1}\over
Q↓0{+}Q↓1z{+}\cdots{+}Q↓{N-1}z↑{N-1}} = W↓N + W↓{N+1}z
+\cdots + W↓{2N-1}z↑{N-1} + O(z↑N)$$
and return to step N2.\quad\blackslug

\yyskip The running time for this
algorithm to obtain the coefficients up to $N = 2↑K$ is $T(N)$,
where
$$T(2N) = T(N) +\hbox{(time to do step N4)} + O(N).\eqno (17)$$
Straightforward algorithms for composition and
division in step N4 will take order $N↑3$ steps, so Algorithm
N will run slower than Algorithm T\null. However, Brent and Kung
have found a way to do the required composition of power series
with $O(N\log N)↑{3/2}$ arithmetic operations, and exercise
6 gives an even faster algorithm for division; hence (17) shows
that power series reversion can be achieved by doing only $O(N
\log N)↑{3/2}$ operations as $N → ∞$.\xskip (On the other hand the
constant of proportionality is such that $N$ must be really
large before Algorithms L and T lose out to this ``high-speed''
method.)

{\sl Historical note:} J. N. Bramhall and M. A.
Chapple published the first $O(n↑3)$ method for power series
reversion in {\sl CACM \bf 4} (1961), 317--318, 503. It was an
\hbox{off-line} algorithm whose running time was approximately the
same as that of Algorithms L and T.

\subsectionbegin{Iteration of series} If we want to study the behavior of an
iterative process $x↓n←f(x↓{n-1})$, we are interested in studying the $n$-fold
composition of a given function $f$ with itself, namely $x↓n=f\biglp f(\ldotsm
f(x↓0)\ldotsm)\bigrp$. Let us define $f↑0(x)=x$ and $f↑n(x)=f\biglp
 f↑{n-1}(x)\bigrp$, so that
$$f↑{m+n}(x)=f↑m\biglp f↑n(x)\bigrp\eqno(18)$$
for all integers $m$, $n≥0$. It also makes sense to define $f↑n(x)$ when $n$
is a negative integer, namely to let $f↑n$ and $f↑{-n}$ be inverse functions
such that $x=f↑n\biglp f↑{-n}(x)\bigrp$; then (18) holds for {\sl all\/} integers
$m$ and $n$. Reversion of series is essentially the operation of finding the
inverse function $f↑{-1}(x)$; for example, Eqs.\ (10) and (11) essentially state
that $z=V\biglp W(z)\bigrp$ and that $t=W\biglp V(t)\bigrp$, so $W=V↑{-1}$.

Suppose we are given two power series $V(z)=z+V↓2z↑2+\cdots$ and $W(z)=z+W↓2z↑2+
\cdots$ such that $W=V↑{-1}$. Let $u$ be any nonzero constant, and consider the
function
$$U(z)=W\biglp u V(z)\bigrp.\eqno(19)$$
It is easy to see that $U\biglp(U(z)\bigrp=W\biglp u↑2V(z)\bigrp$, and in general
that
$$U↑n(z)=W\biglp u↑nV(z)\bigrp\eqno(20)$$
for all integers $n$. Therefore we have a simple expression for the $n$th iterate
$U↑n$, which can be calculated with roughly the same amount of work for all $n$.
Furthermore, we can even use (20) to define $U↑n$ for non-integer values of $n$;
the ``half iterate'' $U↑{1/2}$, for example, is a function such that $U↑{1/2}\biglp
U↑{1/2}(z)\bigrp=U(z)$.\xskip$\biglp$There are two such functions $U↑{1/2}$,
obtained by using $\sqrt u$ and $-\sqrt u$ as the value of $u↑{1/2}$.$\bigrp$

We obtained the simple state of affairs in (20) by starting with $V$ and $u$, then
defining $U$. But in practice we generally want to go the other way: Starting
with some given function $U$, we want to find $V$ and $u$ such that (19) holds,
i.e., such that
$$V\biglp U(z)\bigrp\,=\,u\,V(z).\eqno(21)$$
Such a function $V$ is called the {\sl Schr\"oder function} of $U$, because it
was introduced by Ernst Schr\"oder in {\sl Math.\ Annalen \bf3} (1871), 296--322.
Let us now look at the problem of finding the Schr\"oder function
$V(z)=z+V↓2z↑2+\cdots$ of a given power series $U=U↓1z+U↓2z↑2+\cdotss$. Clearly
$u=U↓1$ if (21) is to hold.

Expanding (21) with $u=U↓1$ and equating coefficients of $z$ leads to a sequence
of equations that begins
$$\baselineskip15pt
\eqalign{U↓1↑2V↓2+U↓2⊗=U↓1V↓2\cr
U↓1↑3V↓3+2U↓1U↓2V↓2+U↓3⊗=U↓1V↓3\cr
U↓1↑4V↓4+3U↓1↑2U↓2V↓3+2U↓1U↓3V↓2+U↓2↑2V↓2+U↓3⊗=U↓1V↓4\cr}$$
and so on. Clearly there is no solution when $U↓1=0$ (unless trivially
$U↓2=U↓3=\cdots=0$); otherwise there is a unique solution unless $U↓1$ is a root
of unity. We might have expected that something funny would happen when $U↓1↑n=1$,
since Eq.\ (20) tells us that $U↑n(z)=z$ if the Schr\"oder function exists. For the
moment let us assume that $U↓1$ is nonzero and not a root of unity; then the
Schr\"oder function does exist, and the next question is how to compute it without
doing too much work.

The following procedure has been suggested by R. P. Brent and J. F. Traub. Equation
(21) leads to subproblems of a similar but more complicated form, so we set
ourselves a more general task whose subtasks have the same form: Let us try to find
$V(z)=V↓0+V↓1z+\cdots+V↓{n-1}z↑{n-1}$ such that
$$V\biglp U(z)\bigrp=W(z)V(z)+S(z)+O(z↑n),\eqno(22)$$
given $U(z)$, $W(z)$, $S(z)$, and $n$, where $n$ is a power of 2 and $U(0)=0$. If
$n=1$ we simply let $V↓0=S(0)\vcenter{\hbox{\:@\char'16}}
\biglp1-W(0)\bigrp$, with $V↓0=1$ if $S(0)=0$ and
$W(0)=1$. Furthermore it is possible to go from $n$ to $2n$: First we find $R(z)$
such that
$$V\biglp U(z)\bigrp=W(z)V(z)+S(z)-z↑nR(z)+O(z↑{2n}).\eqno(23)$$
Then we compute
$$\A W(z)=W(z)\biglp z/U(z)\bigrp↑n,\qquad\A S(z)=R(z)\biglp z/U(z)\bigrp↑n\eqno(24)
$$and find $\A V(z)=V↓n+V↓{n+1}(z)+\cdots+V↓{2n-1}z↑{n-1}$ such that
$$\A V\biglp U(z)\bigrp=\A W(z)\A V(z)+\A S(z)+O(z↑n).\eqno(25)$$
It follows that $V↑\ast\biglp U(z)\bigrp=W(z)V↑\ast(z)+S(z)+O(z↑{2n})$ as desired,
when $V↑\ast(z)$ is equal to $V(z)+z↑n\A V(z)$.

The running time $T(n)$ of this procedure satisfies
$$T(2n)=2T(n)+C(n),\eqno(26)$$
where $C(n)$ is the time to compute $R(z)$, $\A W(z)$, and $\A S(z)$. The function
$C(n)$ is dominated by the time to compute $V\biglp U(z)\bigrp$ modulo $z↑{2n}$,
and $C(n)$ presumably grows faster than order $n↑{1+ε}$; therefore the solution
$T(n)$ to (26) will be of order $C(n)$. For example, if $C(n)=cn↑3$ we have
$T(n)\approx{4\over3}cn↑3$; or if $C(n)$ is $O(n\log n)↑{3/2}$ using ``fast''
composition, we have $T(n)=O(n\log n)↑{3/2}$.

The procedure breaks down when $W(0)=1$ and $S(0)≠0$, so we need to investigate
when this can happen. It is easy to prove by induction on $n$ that the solution of
(22) by the Brent--Traub method entails consideration of exactly $n$ subproblems,
in which the coefficient of $V(z)$ on the right-hand side takes the respective 
values $W(z)\biglp z/U(z)\bigrp\null↑j$ for $0≤j<n$ in some order. If $W(0)=U↓1$ and
if $U↓1$ is not a root of unity, we therefore have $W(0)=1$ only when $j=1$; the
procedure will fail in this case only if (22) has no solution for $n=2$.

The Schr\"oder function for $U$ can therefore be found by solving (22) for
$n=2$, 4, 8, 16, $\ldotss$, with $W(z)=U↓1$ and $S(z)=0$, whenever $U↓1$ is nonzero
and not a root of unity.

If $U↓1=1$, there is no Schr\"oder function unless $U(z)=z$. But Brent
and Traub have found a fast way to compute $U↑n(z)$ even when $U↓1=1$, by
making use of a function $V(z)$ such that
$$V\biglp U(z)\bigrp=U↑\prime(z)V(z).\eqno(27)$$
If two functions $U(z)$ and $\A U(z)$ both satisfy (27), for the same $V$, it
is easy to check that their composition $U\biglp\A U(z)\bigrp$ does too; therefore
all iterates of $U(z)$ are solutions of (27). Suppose that $U(z)=z+U↓{k\,}z↑k+U↓{k+1}
z↑{k+1}+\cdotss$ where $k≥2$ and $U↓k≠0$. Then it can be shown that there is a
unique power series of the form $V(z)=z↑k+V↓{k+1}z↑{k+1}+V↓{k+2}z↑{k+2}+\cdots$
satisfying (27). Conversely if such a function $V(z)$ is given, and if $k≥2$ and
$U↓k$ are given, then there is a unique power series of the form $U(z)=z+U↓{k\,}z↑k
+U↓{k+1}z↑{k+1}+\cdots$ satisfying (27). The desired iterate $U↑n(z)$ is the
unique power series $P(z)$ satisfying
$$V\biglp P(z)\bigrp=P↑\prime(z)V(z)\eqno(28)$$
such that $P(z)=z+nU↓{k\,}z↑k+\cdotss$. Both $V(z)$ and $P(z)$ can be found by
appropriate algorithms (see exercise 13).

If $U↓1$ is a $k$th root of unity, but not equal to 1, the same method can be
applied to the function $U↑k(z)=z+\cdotss$, and $U↑k(z)$ can be found from
$U(z)$ by doing $l(k)$ composition operations (cf.\ Section 4.6.3). We can also
handle the case $U↓1=0$: If $U(z)=U↓{k\,}z↑k+U↓{k+1}z↑{k+1}+\cdots$ where $k≥2$ and
$U↓k≠0$, the idea is to find a solution to the equation
$V\biglp U(z)\bigrp=U↓kV(z)↑k$; then $U↑n(z)=V↑{-1}\biglp U↓k↑{(k↑n-1)/(k-1)}V(z)
↑{k↑n}\bigrp$. Finally, if $U(z)=U↓0+U↓1z+\cdots$ where $U↓0≠0$, let $α$ be a
``fixed point'' such that $U(α)=α$, and let $\A U(z)=U(α+z)-α=zU↑\prime(α)+z↑2
U↑{\prime\prime}(α)/2!+\cdotss$; then $U↑n(z)={\A U}↑n(z-α)+α$. Further details
can be found in Brent and Traub's paper [{\sl SIAM J. Computing}, to appear].

\exbegin{EXERCISES}

\exno 1. [M10] The text explains
how to divide $U(z)$ by $V(z)$ when $V↓0 ≠ 0$; how should the
division be done when $V↓0 = 0?$

\exno 2. [20] If the coefficients of $U(z)$ and $V(z)$ are integers
and $V↓0 ≠ 0$, find a recurrence relation for the integers $V↑{n+1}↓{0}W↓n$,
where $W↓n$ is defined by (3). How would you use this for power
series division?

\exno 3. [M15] Does formula (9) give the right results when
$α = 0$? When $α = 1$?

\trexno 4. [HM23] Show that simple modifications of (9) can be
used to calculate $e↑{V(z)}$ and $\ln\biglp 1 + V(z)\bigrp$,
when $V(z) = V↓1z + V↓2z↑2 +\cdotss$.

\exno 5. [M00] What happens when a power series is reverted
twice; i.e., if the output of Algorithm L or T is reverted again?

\trexno 6. [M21] (H. T. Kung.)\xskip Apply Newton's method to the computation
of $W(z) = 1/V(z)$, when $V(0) ≠ 0$, by finding the power-series
root of the equation $f(x)=0$, where $f(x) = x↑{-1} - V(z)$.

\exno 7. [M23] Use Lagrange's inversion formula (12) to find
a simple expression for the coefficient $W↓n$ in the reversion
of $z = t - t↑m$.

\trexno 8. [M25] Lagrange's inversion formula can be generalized
as follows: If $W(z) = W↓1z + W↓2z↑2 +\cdots = G↓1t
+ G↓2t↑2 + G↓3t↑3 +\cdots = G(t)$, where $z$ and
$t$ are related by Eq.\ (10), then $W↓n = U↓{n-1}/n$ where
$$U↓0 + U↓1z + U↓2z↑2 +\cdots = (G↓1 + 2G↓2z + 3G↓3z↑2
+\cdotss)(1 + V↓2z + V↓3z↑2 +\cdotss)↑{-n}.$$
(Equation (12) is the special case $G↓1 = 1$, $G↓2
= G↓3 =\cdots = 0$. This equation can be proved,
for example, by using tree-enumeration formulas as in exercise
2.3.4.4--33.)
Extend Algorithm L so that it obtains the coefficients
$W↓1$, $W↓2$, $\ldotss$, in this more general situation, without
substantially increasing its running time.

\exno 9. [11] Find the values of $T↓{mn}$
computed by Algorithm T as it determines the first five coefficients
in the reversion of $z = t - t↑2$.

\exno 10. [M20] Given that $y = x↑α + a↓1x↑{α+1}+a↓2x↑{α+2} +\cdotss$,
$α ≠ 0$, show how to compute the coefficients in the expansion
$x = y↑{1/α} + b↓2y↑{2/α} + b↓3y↑{3/α} +\cdotss$.

\trexno 11. [M25] Let
$$U(z) = U↓0 + U↓1z + U↓2z↑2 +\cdots \quad\hbox{and}\quad
V(z) = V↓1z + V↓2z↑2 + V↓3z↑3 +\cdotss;$$
design an algorithm that computes the first $N$
coefficients of $U\biglp V(z)\bigrp$.

\exno 12. [M20] Find a connection between polynomial division
and power-series division: Given polynomials $u(x)$ and $v(x)$
of respective degrees $m$ and $n$ over a field, show how to
find the polynomials $q(x)$, $r(x)$ such that $u(x) = q(x)v(x)
+ r(x)$, deg$(r) < n$, using only operations on power series.

\trexno 13. [HM30] Fill in the details of Brent and Traub's method for
calculating $U↑n(z)$ when $U(z)=z+U↓{k\,}z↑k+\cdotss$, using (27) and (28).

\vfill\eject % This makes the following quotation appear on a right-hand page
\vfill
\quoteformat{And it shall be,\cr
when thou hast made an end of reading this book,\cr
that thou shalt bind a stone to it,\cr
and cast it into the midst of Euphrates.\cr}
\author{Jeremiah 51:63}
\eject